package com.github.zyf.ioaj.algorithm.physics;

import com.github.zyf.ioaj.algorithm.AbstractAlgorithm;
import com.github.zyf.ioaj.algorithm.Input;
import com.github.zyf.ioaj.algorithm.Output;
import com.github.zyf.ioaj.annotation.internal.Algorithm;
import com.github.zyf.ioaj.annotation.internal.BasicParam;
import com.github.zyf.ioaj.annotation.internal.TempParam;
import com.github.zyf.ioaj.constant.FileType;
import com.github.zyf.ioaj.constant.ParamType;
import com.github.zyf.ioaj.export.Report;
import com.github.zyf.ioaj.export.Tracker;
import com.github.zyf.ioaj.utils.basedata.ArrayUtils;

/**
 * <h1>Gravitational Search Algorithm</h1>
 * <p>引力搜索算法</p>
 * <p>默认求解最小值，如果需要求解最大值，将目标函数乘以-1即可，最后得到的值再次乘以-1就是得到的最优值（自变量取值不用乘以-1）</p>
 * <p>需要计算则调用{@link #go()}</p>
 * <p>需要跟踪每次迭代则调用{@link #tracker(String, String, String, FileType)}</p>
 * <p>需要生成报告则调用{@link #report(String, String, String, FileType)}</p>
 * <p>这种算法的灵感来自于牛顿的万有引力定律与运动定律：</p>
 * <ol>
 *     <li>任意两个质点有通过连心线方向上的力相互吸引，该引力大小与它们质量的乘积成正比与它们距离的平方成反比。</li>
 *     <li>力使物体获得加速度。</li>
 * </ol>
 * <p>
 * 在GSA中，质点被抽象成解空间的一个解，解之间存在一个相互吸引力，这个吸引力由解的质量与两个解之间的距离确定，
 * 质点的质量被抽象成解的评估函数值。在解空间中，每个解由其他解对其的吸引力获得加速度，质量更大（评估函数值更优）所提供的加速度更大，
 * 从而使解向更优解的方向移动。
 * </p>
 * <p>由此，GSA和PSO算法两者有很多相似之处，但是GSA与PSO的外层框架是一致的，它们最关键的粒子移动策略不同：</p>
 * <ol>
 *     <li>在PSO中，粒子的移动只使用两个最佳位置pbest和gbest。但是在GSA中，智能体的方向是根据所有其他智能体或部分较优的智能体获得的总力来计算的。</li>
 *     <li>PSO使用一种内存来更新速度（由于pbest和gbest）。然而，GSA是无内存的，只有智能体的当前位置在更新过程中起作用。</li>
 *     <li>在PSO中，更新不考虑解之间的距离，而在GSA中，力与解之间的距离成反比。</li>
 *     <li>PSO模拟了鸟类的社会行为，而GSA的灵感来自于物理现象。</li>
 * </ol>
 * <p>由此可以概括基本步骤</p>
 * <ol>
 *     <li>确定搜索空间。</li>
 *     <li>随机初始化每个粒子。</li>
 *     <li>对粒子进行适应度评价。</li>
 *     <li>更新引力常量G，更新粒子中最好的个体 best与最差的个体worst，更新粒子的质量。</li>
 *     <li>计算每个粒子在不同方向上的总引力。</li>
 *     <li>计算粒子的加速度，基于此更新粒子的速度。</li>
 *     <li>根据速度更新粒子在解空间中的位置。</li>
 *     <li>重复步骤3-7直到达到停止标准。</li>
 * </ol>
 *
 * @author zyf
 */
@Algorithm
public class GSA extends AbstractAlgorithm {
    /**
     * 粒子个数
     */
    @BasicParam(type = ParamType.INT, desc = "粒子个数")
    private int n;

    /**
     * 粒子的加速度
     */
    @TempParam(type = ParamType.DOUBLE_ARRAY_2, desc = "粒子的加速度")
    private double[][] acceleration;

    /**
     * 粒子的位置
     */
    @TempParam(type = ParamType.DOUBLE_ARRAY_2, desc = "粒子的位置")
    private double[][] position;

    /**
     * 粒子的速度
     */
    @TempParam(type = ParamType.DOUBLE_ARRAY_2, desc = "粒子的速度")
    private double[][] velocity;

    /**
     * 粒子的质量
     */
    @TempParam(type = ParamType.DOUBLE_ARRAY_1, desc = "粒子的质量")
    private double[] mass;

    /**
     * 粒子的质量比
     */
    @TempParam(type = ParamType.DOUBLE_ARRAY_1, desc = "粒子的质量比")
    private double[] Mass;

    /**
     * 每个粒子的适应度
     */
    @TempParam(type = ParamType.DOUBLE_ARRAY_1, desc = "粒子适应度")
    private double[] fitness;

    /**
     * 引力常量初值
     */
    @BasicParam(type = ParamType.DOUBLE, desc = "引力常量初值")
    private double G0;

    /**
     * 引力常量
     */
    @TempParam(type = ParamType.DOUBLE, desc = "引力常量")
    private double G;

    /**
     * 当前最优解
     */
    @TempParam(type = ParamType.DOUBLE_ARRAY_1, desc = "当前最优解")
    private double[] x;

    /**
     * 抽象算法构造函数，需要提供输入量，即提前初始化好{@link Input}对象
     *
     * @param input {@link Input}
     */
    public GSA(Input input) {
        super(input);
        init();
    }

    private void init() {
        this.n = getIter();
        this.G0 = 92;
        this.G = G0;

        this.x = new double[getInput().getVarNum()];
        this.acceleration = new double[n][getInput().getVarNum()];
        this.position = new double[n][getInput().getVarNum()];
        this.velocity = new double[n][getInput().getVarNum()];

        this.mass = new double[n];
        this.Mass = new double[n];
        this.fitness = new double[n];

        for (int i = 0; i < n; i++) {
            for (int j = 0; j < getInput().getVarNum(); j++) {
                position[i][j] = getInput().randomValueInRange(j);
                velocity[i][j] = 0.2 * getInput().randomValueInRange(j);
            }
        }
        x = position[0].clone();
        calculateMass();
        calculateAcceleration();
    }

    /**
     * 每次迭代时扫描一次
     *
     * @param t 迭代次数
     */
    public void scan(int t) {
        G = G * Math.exp(-20.0 * t / getIter());
        calculateMass();
        calculateAcceleration();
        updateVelocityAndPosition();
    }

    /**
     * 生成最优解
     */
    private void generateOptimization() {
        // 获取最优值
        double optimalOutput = getInput().getFunction().apply(x);
        setOutput(new Output(optimalOutput, x));
    }

    /**
     * 运行并计算结果
     */
    @Override
    public void go() {
        int t;
        for (t = 0; t < getIter(); t++) {
            scan(t);
        }
        generateOptimization();
    }

    /**
     * 运行并计算结果（按迭代逐步跟踪，自带运行，无需手动），如果目的是计算出结果，请使用{@link #go()}
     *
     * @param author   作者
     * @param message  备注
     * @param title    标题
     * @param fileType 文件类型 {@link FileType}
     * @return {@link Tracker}
     */
    @Override
    public Tracker tracker(String author, String message, String title, FileType fileType) {
        int t;
        Tracker tracker = new Tracker(this, author, message, title, fileType);
        tracker.startTracker();
        for (t = 0; t < getIter(); t++) {
            scan(t);
            tracker.trackOne(t + 1);
        }
        generateOptimization();
        tracker.endTracker();
        return tracker;
    }

    /**
     * 生成报告（要先手动运行才能生成报告）
     *
     * @param author   作者
     * @param message  备注
     * @param title    标题
     * @param fileType 文件类型 {@link FileType}
     * @return {@link Report}
     */
    @Override
    public Report report(String author, String message, String title, FileType fileType) {
        return new Report(this, author, message, title, fileType);
    }

    /**
     * 根据适应度计算质量比（类似于归一化）
     */
    private void calculateMass() {
        int i;
        // 先找到最好和最差的
        double max = ArrayUtils.max(fitness);
        double min = ArrayUtils.min(fitness);
        // 如果最大值和最小值一样了，所有的质量就都设置为1
        if (max == min) {
            for (i = 0; i < n; i++) {
                mass[i] = 1.0;
            }
        } else {
            for (i = 0; i < n; i++) {
                mass[i] = (fitness[i] - max) / (min - max);
            }
        }
        // 计算质量比
        double sum = ArrayUtils.sum(mass);
        for (i = 0; i < mass.length; i++) {
            Mass[i] = mass[i] / sum;
        }
    }

    public void calculateAcceleration() {
        double[] Fi = new double[n];

        int i, j, k;
        for (i = 0; i < n; i++) {
            for (j = 0; j < n; j++) {
                if (i != j) {
                    for (k = 0; k < getInput().getVarNum(); k++) {
                        Fi[i] += Math.random() * G * ((Mass[i] * Mass[j])
                                / (calcDistance(position[i], position[j]) + 1))
                                * (position[j][k] - position[i][k]);
                    }
                }
            }
            if (Mass[i] != 0) {
                for (k = 0; k < getInput().getVarNum(); k++) {
                    acceleration[i][k] = Fi[i] / Mass[i] / 10;
                }
            } else {
                for (k = 0; k < getInput().getVarNum(); k++) {
                    acceleration[i][k] = 10.0;
                }
            }
            for (k = 0; k < getInput().getVarNum(); k++) {
                Fi[i] = 0.0;
            }
        }
    }

    /**
     * 计算速度和位置
     */
    private void updateVelocityAndPosition() {
        int i, j;
        for (i = 0; i < n; i++) {
            for (j = 0; j < getInput().getVarNum(); j++) {
                velocity[i][j] = Math.random() * velocity[i][j] + acceleration[i][j];
                position[i][j] = position[i][j] + velocity[i][j];

                // 防止越界
                if (position[i][j] < getInput().getLb()[j]) position[i][j] = getInput().getLb()[j];
                if (position[i][j] > getInput().getUb()[j]) position[i][j] = getInput().getUb()[j];
            }
            fitness[i] = getInput().getFunction().apply(position[i]);
            if (fitness[i] < getInput().getFunction().apply(x)) {
                x = position[i].clone();
            }
        }
    }

    /**
     * 计算距离（欧几里得距离）
     *
     * @param p1 点1
     * @param p2 点2
     * @return double
     */
    private double calcDistance(double[] p1, double[] p2) {
        int j;
        double d2 = 0.0;
        for (j = 0; j < getInput().getVarNum(); j++) {
            d2 += (p1[j] - p2[j]) * (p1[j] - p2[j]);
        }
        return Math.sqrt(d2);
    }

    public int getN() {
        return n;
    }

    public void setN(int n) {
        this.n = n;
    }

    public double getG0() {
        return G0;
    }

    public void setG0(double g0) {
        G0 = g0;
    }

    public double[][] getAcceleration() {
        return acceleration;
    }

    public double[][] getPosition() {
        return position;
    }

    public double[][] getVelocity() {
        return velocity;
    }

    public double[] getMass() {
        return mass;
    }

    public double[] getFitness() {
        return fitness;
    }

    public double getG() {
        return G;
    }

    public double[] getX() {
        return x;
    }
}
